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Abstract 

We determine for QCD with nj = 2 dynamical quark flavors by fitting the QQ static 
potential known analytically in the perturbative regime up to terms of 0{ag) and ~ q^Iuq^ 
to corresponding results obtained from lattice simulations. This has become possible, due to 
recent advances in both perturbative calculations, namely the determination and publication of 
the last missing contribution to the QQ static potential at ©(a^), and lattice simulations with 
rif = 2 dynamical quark flavors performed at the rather fine lattice spacing of a ~ 0.042 fm. 
Imposing conservative error estimates we obtain Aj;^ = 315(30) MeV. 



1 Introduction 



The last 15 years have seen substantial progress in the determination of the QQ (quark- 
antiquark) static potential in the perturbative regime. Whereas it had been known up to 
0{a1) for quite a long time [Il[2], its next order contribution ~ was determined only about 
15 years ago, in the end of the 1990s [31 HI [5]. It took another decade until all contributions at 
0(af ) were available in the MS-scheme [3 EJ O [TU]. A strict power-series expansion in Ug, 
however, is not the whole story. Notably, already in the early days of quantum chromodynamics 
(QCD), it had been recognized jTT] that the static potential does not have a strict power-series 
expansion in a^. Beyond 0(ag), also logarithmic contributions in are induced. The leading 
one ~ In has been determined explicitly, and in a systematic way by [12] . In summary, 
the static potential is now completely known up to C'(af) and ~ Ina^- On the other hand, 
the QQ static potential is easily accessible in lattice QCD simulations (cf. e.g. [I3l[l3]), and has 
received continuous attention since the early days of lattice QCD. Having focused on the study 
of pure Yang-Mills theory without dynamical quark flavors, i.e. n/ = 0, in the beginning, it has 
meanwhile become common to also account for light-flavor degrees of freedom. Moreover, the 
increase in available computing power has opened up the possibility to simulate on larger and 
larger lattices, thereby allowing for smaller and smaller lattice spacings. 

Whereas the aforementioned perturbative calculations are conventionally performed in momen- 
tum space, lattice QCD naturally provides the static potential in position space. Due to the fact 
that QCD is asymptotically free, perturbation theorjH is viable at large momentum-transfers. 
After Fourier transform to position space, this allows for trustworthy insights at small QQ sep- 
arations. In contrast, lattice simulations at a given lattice spacing cannot resolve arbitrarily 
small separations. The refinements on the lattice as well as the perturbative side have signifi- 
cantly reduced, and recently even bridged, the gap between the ranges of applicability of both 
approaches. This allows for a quantitative comparison of perturbative calculations with lattice 
simulations at QQ separations, where both approaches are expected to be applicable. While first 
exploratory attempts date back to the early 1990s [15], recently a determination [16] of ?^o^ms 
= has been performed, resorting to the effective field theory of potential non-relativistic 
QCD (pNRQCD) [T5] and the lattice results of [13]. In this paper we focus on QCD with 
Uf = 2 dynamical quark flavors and aim at the determination of Aj^ from the static potential. 
In contrast to the aforementioned study at ny- = 0, we do not explicitly resort to pNRQCD, but 
directly confront the perturbatively determined QQ static potential up to terms of 0{a1) and 
~ In Us with lattice results. 

Our paper is organized as follows. Whereas section [2] outlines the determination of the QQ 
static potential by means of lattice QCD, section [3] focuses on its evaluation within perturbation 
theory. In section U] we describe, how Aj^ and the associated systematic and statistical errors 
are determined by fitting various orders of the perturbative expansion to lattice results of the QQ 
static potential. We briefly compare our result for Aj^ at n/- = 2 to previous determinations, 
and end with conclusions in section [5l 



In agreement with the prevalent notation, particularly in the field of lattice QCD, we use the terms static 
potential and static energy synonymously. Note, however, that sometimes a distinction is made and these terms 
refer to different quantities. 

^We use the term perturbation theory in the sense, that the coupling as is small and can be dealt with 
perturbatively. 
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2 The static potential from lattice QCD 



Although we compute the static potential for QCD with Uf = 2 dynamical quark flavors, 
employing standard methods from lattice gauge theory, we will give a description of our lattice 
calculation below to make the paper self-contained. 

2.1 Lattice setup 

We use Uf = 2 gauge link configurations generated by the European Twisted Mass Collaboration 
(ETMC) [IHl dni [21] . The gauge action is tree- level Symanzik improved [22] . 



Sg[U] = 




with 6o = 1 — 86i and hi = —1/12. The quark action is Wilson twisted mass (cf. [23 ] [24 ] [25 ] [26] ) . 
Sf[x,X,U] = a^^x{x){Dw + ilJ'(ilbTz)x{x), (2) 

with 

= ]^ (7^ (V^ + V;) - aV* V^) + mo. (3) 

Here a denotes the lattice spacing, and V* are the gauge covariant forward and backward 
derivatives, mo and /iq are the bare untwisted and twisted quark masses, respectively, rs is the 
third Pauli matrix acting in flavor space, and x = ix^^\ X^'^'') represents the quark flelds in the 
so-called twisted basis. The twist angle uj is given by tanu = /in/mR, where /xr and mR denote 
the renormalized twisted and untwisted quark masses, uj has been tuned to tt/2 by adjusting 
mo appropriately (cf. [20] for details). This ensures automatic 0{a) improvement for many 
observables, including the static potential. 

The ensembles of gauge link configurations considered here are collected in Table [TJ We have 
data for four different values of the lattice spacing, a wide range of pion masses mps and 
spacetime volumes x T. Details on the generation of these gauge field configurations as well 
as on the computation and the analysis of standard quantities (e.g. lattice spacing or pion mass) 
can be found in [201 IZI]- We also provide the number of gauge link configurations, used in the 
computation of the static potential, for each ensemble. 

2.2 Computation of the static potential 

As usual, the static potential V{r) at QQ separation r = |r| is extracted from the exponential 
decay of Wilson loop averages {W{r, t)) with respect to their temporal extent t, while keeping 
their spatial extent r constant. We consider Wilson loops formed by APE smeared spatial links 
(A^APE = 60, ttAPE = 0.5 for all our gauge link ensembles; cf. for details), and ordinary, i.e. 
unsmeared temporal links. 
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P 


r 

a m im 


[j-'/ci) X 1 /a 


mps in MeV 


# gauges 


3.90 


0.079(3) 


24^ X 48 


340(13) 


168 


4.05 


0.063(2) 


323 X 64 


325(10) 


71 








449(14) 


100 








517(16) 


92 


4.20 


0.0514(8) 


24^ X 48 


284(5) 


123 






48^ X 96 




46 


4.35 


0.0420(17) 


32^ X 64 


352(22) 


146 



Table 1: Ensembles of gauge link configurations; the scale has been set via the pion mass and 
the pion decay constant, using chiral perturbation theory (cf. e.g. [201 121j). 



Given nowadays typical lattice spacings of >0.04fm, there is only a small range of QQ separa- 
tions, where the static potential V{r) can be computed reliably by means of both lattice QCD 
and perturbation theory. For our finest lattice spacing (/3 = 4.35 ensemble) this range roughly 
amounts to 3o ... 5a (cf. section HT3l for a detailed discussion). Therefore, it is desirable to obtain 
a large number of lattice data points in this interval, i.e. a resolution significantly finer than 
the lattice spacing. To this end, we do not only consider ordinary on-axis, but also so-called 
off- axis Wilson loops, whose spatial sides are not exclusively oriented parallel to the x-, y- or 
z-axis, respectively. The off-axis sides are constructed by optimally approximating a straight 
line through a product of on-axis links and/or two- and three-dimensional diagonal links. In 
terms of on-axis links, two dimensional diagonal links are defined as 

Uix;x + fi + i)) = Psu{3){u^{x)U^{x + il) + U^{x)Uf,{x + i))y (4) 
with the projector -Pgj^jg), implementing a projection to SU{3), 

Psuism = , U' = U{Wuf\ (5) 

Three dimensional diagonal links are computed analogously, 

U{x;x + fi + u + p) = 

= Psu{z)(yi,{x)Uy{x + jX)Up{x + iX + v) + U^{x)Up{x + fi)U^{x + fi + p) + ...y (6) 

The optimal approximation of the straight line in terms of the above mentioned links is con- 
structed efficiently by employing the Bresenham algorithm [28], which has already been applied 
to the computation of the QQ static potential in previous works (cf. e.g. [29]). Following this 
approach, we have computed all Wilson loop averages (VF(|r|,i)) with r = na, n € and 
\r\ < 10a. To improve the signal quality, we have averaged over loops, which are related by 
discrete translational and cubic rotational symmetry. 

The QQ static potential V{r) at separation r is obtained in a two-step procedure. First, an 
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effective potential is computed, 



T.(effective). .^ _ 1 / {W{r,t)) \ 



(7) 



Second, the t-independent quantity V{r) is obtained by performing an uncorrelated minimiz- 
ing fit to y(effective)(^^^) 

in a suitable t range. This range is chosen such that excited states are 
strongly suppressed, while statistical errors are still small. 

Even though physical observables are automatically 0{a) improved, as we are employing Wilson 
twisted mass lattice QCD at maximal twist (cf. section 12. ip . we reduce lattice discretization 
errors even further. To this end, we use a method to improve the static potential, which is 
explained in [13l[30]. The improved potential is defined as follows, 

y (improved) i^^(improvcd) -J _ (g) 

with 7- (improved) determined from 



4^^(improved) (2vr)3 J _^ 4 ^3^^ sin2(fej/2) + (4/3) S^Li sm\kj/2) 

(cf. [30], section 3.2.2). This procedure assures that at tree level the lattice static potential 
computed with the tree-level Symanzik improved action is identical to the continuum static 
potential, i.e. y(ii^P™ved) 

is a tree-level improved observable. \/(i™pi'oved) -g ^j^^ quantity to be 
confronted with perturbative calculations in section [H 

We have solved the integral ([9]) numerically via standard Monte Carlo sampling. Decomposing 
the integral and treating its singular part analytically, we were able to reach a precision, where 
statistical errors are negligible (cf. appendix!^ for details). 

Figure [J compares the unimproved (left plot) and the improved (right plot) static potential for 
P = 4.35. In the unimproved static potential, discretization errors are clearly visible at small 
QQ separations, where they amount to a strong violation of rotational invariance. On the other 
hand, the improved potential is rather smooth and no such violation is visible. 
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unimproved static potential (/3 = 4.35) 



improved static potential (/3 = 4.35) 




1 1.5 



2 2.5 3 3.5 

r/a 



1 1.5 



2 2.5 3 3.5 4 

r/a 



Figure 1: Comparison of the unimproved (left) and the improved (right) static potential for 
/3 = 4.35; the solid straight lines connecting the data points are drawn to guide the eye. 
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3 The static potential in perturbation theory 



We now turn to the QQ static potential in perturbation theory. It is completely known up 
to terms of order and the contribution ~ (/z) In where ^ refers to an a priori 

arbitrary, but fixed momentum scale in the perturbative regime, i.e. /i ^ Aqcd- As we are 
interested in QCD, we limit our considerations to the color gauge group SU(3). Instead of 
referring to the corresponding literature only, we aim at keeping this paper rather self-contained. 
Therefore, we explain and explicitly provide all the formulae relevant for this work here. 

As a matter of fact, most of the calculations in perturbative QCD are conveniently performed 
in momentum space. This is also true for the static potential. However, aiming at a comparison 
with lattice simulations naturally performed in position space, we are ultimately interested in a 
perturbative expression of the static potential in position space. 

3.1 The static potential in momentum space 
Defining [31] 



L = L{fi,p) = In 4, (10) 

with p = \p\, the static potential in momentum space can be written as 
47r 

V{p) = -CF^av[as{fi),L{fi,p)], (11) 



where 

ay[as(/i),L(/i,p)] = Usifi) 



+ 



P-iiL) + a3inln asilJ.) +...}. (12) 



All the terms written explicitly in eq. I\12\\ are known. For SU(3) Cp = 4/3. The expansion 
coefficients Pn{L) are polynomials in L [31J and are stated below (cf. eqs. (|20|) to (|22p ). 

Up to 0{a^g{^)) the static potential has a strict power-series expansion in as(/x). Beyond this 
order the power-series expansion in as{n) breaks down [11] and one also encounters logarithmic 
contributions in as{fJ.). The first such term is ~ ag(/i) In as (^). The contribution linear in Us in 
eq. (|12p is the leading order (LO) expression of the static potential. Terms up to 0{a1) corre- 
spond to next-to-leading order (NLO), and up to 0{ag) to next-to-next-leading order (NNLO). 
Finally, terms up to ©(af) together with the Inag contribution are referred to as NNNLO in 
the following. 

Let us emphasize that, whereas we wrote the static potential by means of a perturbative ex- 
pansion in terms of the coupling as(^) at a fixed momentum-scale fi in eqs. (jlip and (jl2p . the 
full static potential should not depend on any externally set scale. It should rather form a 
renormalization group (RG) invariant. This is true order by order in a perturbative expansion 
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in as also. Hence, av[ois, L], exhibiting an implicit /^-dependence via both its arguments as and 
L, has to obey the following renormalization group equation [8], 

d I d a d \ 

IJ'-^av[as{fJ'),L{^i,p)] = ^ ( aL T^^""" J "^f""" ^ ^' ^^^^ 

Here /3[acj(/i)] denotes the QCD /3-function, characterizing the running of the coupling as, defined 
as 

- -J^^. (14) 

It has the following power-series expansion in as{f^), 

n=0 ^ ^ 

with the expansion coefficients Pn known up to n = 3, i.e. to 4-loop order. Whereas /3o and f3i 
are independent of the renormalization scheme, P3 and /34 are scheme-dependent. They have 
been determined for arbitrary compact semi-simple Lie groups in the MS-scheme |32j . For SU(3) 
they read 

/3o = 11-^"/. (16) 
38 

/3i = 102 --nj, (17) 



3 

/32 = rifH nl (18) 

^ 2 18 54 ^ ^ 

/ 149753 / 1078361 6508 / 50065 6472 2 

+ ^4 (19) 

729 ■'^ ^ ' 

with ?ij denoting the number of massless, dynamical quark flavors. 

Eq. (jl3p constrains the polynomials Pn{J->) hi eq. (jl2p to be of the following form [8], 

Pi(L) = ai+/3oL, (20) 
^'2(^) = a2 + (2ai/3o + /3i)L + /32L2, (21) 

P3(^) = a3 + (3a2/3o + 2ai/3i + /32)L + /3o (3ai/3o + ^/3i^ L2 + . (22) 

Prefactors of ag'''^(//) contain terms up to order (n € N). Hence, in the perturbative limit 
considered here the expansion is not only in asi^li), but also in powers of L. For eq. (I12p to 
yield trustworthy results it is, therefore, not enough to require that as{f^) is small. Besides 
fi,p ^ AqcD; in general one has also to ensure that L ~ 1, i.e. that p and fi are of the same 
order. 



2857 5033 325 2 
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For fixed nj, an {n = 1,2,3) and aain are constants. Their values cannot be extracted from 
eq. (fT2]) . but necessitate an explicit determination. Ttie constants oi [H [2] and 02 O [H [5] are 
known analytically. For gauge group SU(3) and in the MS-schema, they read 

31 10 

ai = y-y^/' (23) 



0-2 



/4343 284 A /1229 52 ^, \ 100 n , 



The logarithmic contribution in eq. (jl2p is conveniently determined in position space \12\ [33] . 
Within the effective field theory of potential non-relativistic QCD (pNRQCD) |17[I18| it naturally 
arises by taking into account the ultrasoft (US) contribution to the static potential. However, it 
can also be determined more directly by resumming a certain class of diagrams (Coulomb ladder 
diagrams) in perturbation theory without resorting to pNRQCD |1H I31j. As it also contributes 
to 03, we postpone the specification of the explicit expressions for 03 and a3in, until having 
performed the transition to position space. 



3.2 The static potential in position space 

The static potential in position space is determined by means of a Fourier transform, 

Vir) = J ^e^^"^^F(p). (25) 

In section [3Tl we argued that, for the perturbative expansion of the potential to yield trustworthy 
results in momentum-space, p and /i should be of the same order. Conversely, in eq. (j25p the 
integral is over the full momentum regime. This seems to be a contradiction at first glance 
only. Let us first motivate that an explicit restriction to the perturbative momentum regime, 
p > Aqcd, would basically reproduce eq. (i25]l up to a r-independent global shift term Kjhift, as 

= Iw? '"'^^''^ ~lw? "^''^^ ^^'^^'^ " '"^^ ' ^''^ 

and, because in the perturbative regime Aqcd ^ f , 

" Jw¥ "'"^^"^ 7 ®<''«™ - 1'^* 

As the lattice potential is known up to a global shift only, for our purposes the explicit expression 
for Vshift is not relevant. Second, as the coefficient of a given power of a"(/i), n G N, encodes 
the full momentum dependence at this order, the Fourier transform can be performed order by 
order in as{fi), irrespective of the relative relation between p and /U. By inspection of its result, 
cf. eq. (I29p . one finds that to prevent the respective logarithms from becoming large, one now 
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has to require that 1/r does not deviate much from fi. Hence, the resulting potential in position 
space V{r) can be considered as trustworthy in the regime, where ^ ^ Aqcdj and moreover 
/i and 1/r are of the same order. 



(28) 



Introducing 
L' = L'{fi,r) = ln(A') + 27i?, 

where "Je denotes the Euler-Mascheroni constant, eq. (|25p immediately results in 



V{r 



r Att 



as(/i) 



An 



P2{L') 



+ 



Att 



3 ^ 



+ 



with [31 

Pi{L') 
P2{L') 

PsiL') 



ai + (3oL', 
a2 + (2ai/3o + /3i)L' + 0^ ( L''' + 



IT 



as + (3a2/3o + 2ai/3i + /Ss)^' + ( Sa^fS^ + -/3o/3i ) ( L'^ + — 
+ f3l + vr^L' + 16C(3)] . 



(29) 

(30) 
(31) 

(32) 



Finally, we explicitly specify the coefficients and aain- A comparison with eq. (21) of [3T 
yields 



as + flsin ln(as(/i)) 



16 o^s 
as + y ^ C"! 



In (CAas(/i)) +7£; 



(33) 



where, for SU(3), Ca = 3 and [M] 



as = 
with 



a 



(0) , .(1) 



(2)^2 , ^(3)^S 



a 



(1) 



97 

27ci + — C2, 
lb 

9 5 68993 16624 

— H C4 1 

2 96 81 27 



C(3) + i?C(5). 



(2) _ 93631 ^ 16^4 412 



,(3) 



972 45 
1000 



C(3), 



729 



(34) 

(35) 
(36) 
(37) 
(38) 



The coefficients q {i = 1 ... 4) are only known numerically, ci and C2 have been determined 
independently both in [B] and in [7]. We use the numerical values from [H], who provide smaller 
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statistical error^, 
ci = 502.24(1) , C2 = 
C3 and C4 are given by [TU] . 
C3 = -709.717 , C4 = 
From eq. ()33p it follows that 

16 



-136.39(12). 
-56.83(1). 



03 



In (Ca) + IE 



6 



03 + 1447r^ 



ln3 + 7_B 



16 



TT^C^ = 1447r2. 



(39) 
(40) 

(41) 
(42) 



Therewith, all the coefficients appearing in the perturbative expansion of the static potential, 
eq. ([12]) and ([29|) . respectively, have been assembled. Only the scale /x and in particular the 
functional dependence of as{^i) on ^ have not yet been specified. 

In particular note that eqs. (|14p and (jl5p imply that as(/x) can be expressed as a function of 
as{v) and ln(/x^/z^^). Making a power-series ansatz for q:s(//) and requiring this expression to 
fulfill eqs. ()14p and (jl5p one obtains 



47r 



+ 



Qs(t^) 
47r 



a.s(t^) 
47r 



/5o^l-'$-?/?o/3iln^ + /32|ln^ 



1/2 



1/2 



+ 0(a5(z/)), (43) 



i.e. prefactors of a^''~^(z/) (n € N) contain terms up to order In^ {fi'^ / u'^) . The expansion in 
eq. (j43p is valid in the perturbative regime, i.e. for momenta /i and i/ such that as(/^) and as(i/) 
are small. Moreover, /i and ;/ have to be of the same order, not to spoil the expansion by inducing 
large logarithms. Inserting eq. (j43p in the expressions for the static potential in momentum and 
position space, respectively, and keeping terms up to ©(af) and ~ (z/) lnas(i/), we exactly 
recover eqs. (jlip and (j29p with fj, substituted by u. This clearly confirms that in principle, i.e. 
from a purely theoretical viewpoint detached from any phenomenological considerations, the 
static potential does not explicitly depend on any fixed scale ^u, order by order in a perturbative 
expansion, as demanded and discussed in the context of eq. (fT3]) above. Truncating at a given 
order, nothing favors a fixed given scale as compared to another fixed scale v, fulfilling the 
above requirements. 



3.3 Towards phenomenological applications 

In this section our aim is to make contact with phenomenology. First we line out, how the 
coupling as(//) at the scale n is conventionally fixed within the framework of perturbation 
theory, namely by introducing a generic momentum scale A. Then we discuss, how to choose // 
appropriately. 

■^The errors associated with ci to C4 turn out to be negligible in the context of our determination; therefore, 
we will not discuss them any further. 
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3.3.1 The scale and its relation to the coupling as{fi) 

Recall, that in order to obtain eq. ()43p . we solved eqs. ()14p and (jlSp by invoking a power-series 
ansatz. The result only enabled us to express the coupling at the scale /j, in terms of the coupling 
at the scale v, without allowing us to identify any natural reference scale. Contrarily, the explicit 
integration of eq. (|14p naturally results in the introduction of a reference scale A. Implementing 
the conventional normalization \3A\ [35] and keeping terms beyond one- loop level in eqs. (jl4p and 
(jlSp . the renormalization group invariant parameter A in our conventions reads (cf. e.g. 

1 27r 



A = /i — exp ' ' 



4vr J (3oas{fJ') Jo a. 



+ -L 

P{as) Poas 2/32 



(44) 



If all coefficients beyond /3o in eq. ()15p are neglected, i.e. at one-loop accuracy. 



A = /xexp <j --T—i-T } ■ (45) 



Eqs. (|44p and (j45p allow for the evaluation of as{^j) as a function of the dimensionless ratio /u/A. 
Hence fixing the scale A in principle amounts to knowing as{^) at any momentum scale in the 
perturbative regime. As the QCD /5- function is known at four- loop order only (cf. eq. (|15p ). at 
best as(/i) can be related to /i/A at four-loop accuracy. To obtain as(/x) as a function of ///A at 
two, three and four-loop accuracy, which we denote by a" '°°''(//), with n = 2,3,4 referring to 
the loop accuracy, we insert the respective expressions for /3(a) and solve eq. numerically. 
Only for as~^°°^(^) a simple exact analytical solution can be inferred from eq. (I45p . The quantity 
A is renormalization scheme-dependent. With the expansion coefficients from eqs. (jlSp and (jl9p 
evaluated in the MS-scheme, A = Aj^. 



3.3.2 Setting the scale ^ 



We have manifestly required in eq. (jl3p . and subsequently also shown in the context of eq. (j43p . 
that the static potential does not explicitly depend on any particular choice of the externally 
set scale /u order by order in ag. Hence, in the full all-order expression" of the static potential 
the choice of the scale 11 would be completely arbitrary and the result would not depend on this 
scale at all. 

In this work, however, our aim is to confront the truncated expression, eq. (j29p . with lattice 
results within a certain range of QQ separations, r^\^ < t < rmax) where perturbation theory 
is assumed to be valid. Whenever truncating the series and only taking into account a finite 
number of terms, results certainly depend on the particular choice of the scale In a range, 
where perturbation theory is valid and the requirement to be of the same order as 1/r is fulfilled 
for two different scales and u (cf. the discussion above eq. ([28|) ). the dependence of eq. ([29|) on 
/i and V, respectively, is expected to be rather weak. Moreover, it should become smaller, when 
increasing the order in the perturbative expansion, thereby signalizing a reasonable convergence 
of the perturbative approximation. We demonstrate in section |4] that this is indeed the case. 

Following this philosophy, we fix the scale /i at a value, which is of the same order as 1/r in 
the spatial interval r-^ain ^ r < rmax) where the comparison with lattice results is going to 
be performed. A convenient choice is /x = 2/(r^[^ + r^^.^), which will also be adopted in the 
following section. 
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4 Determination of A 



In this section, we determine in units of the lattice spacing, i.e. aAj^, by fitting perturbative 
expressions for the static potential (cf . section [3|) to corresponding lattice results (cf . section [2]) . 
Using the values of the lattice spacing listed in Table [H these results can easily be converted to 
physical units, e.g. MeV. Unless explicitly specified otherwise, the errors on do not account 
for the uncertainties associated with the lattice spacing. We include these errors at the end, 
when quoting our final result for Aj^ (cf. eq. (08])). 

4.1 Fitting procedures 

Let us for the moment assume that we have already specified the scale /x. To determine Aj^ we 
then apply two strategies: 

(A) First we determine by means of fitting the perturbative expression for the QQ static 
potential in eq. (j29p to the lattice data points. To test the convergence of the perturbative 
expression, and to obtain an estimate for the systematic error, we do this for the four 
different orders, i.e. LO, NLO, NNLO and NNNLO, introduced in section [HTTl In a second 
step we always identify q;s(/u) = at~^°°^(/.i), i.e. employ the QCD /3-function to the best 
accuracy available, and solve eq. ()44p for A^^ (cf. also section r3.3.ip . 

(B) Instead of always identifying as(A*) with as~'°°''(/u) in the perturbative expression for the 
static potential, we explicitly make use of the various expressions relating as(/i) to A^//i 
at different loop accuracies, i.e. ai~^°°''(/i) with I = 1, ...,4. The identification of a<j(A') with 
ai~'°°^(/i) is varied, depending both on the power n of a given contribution ~ a^(//) to the 
static potential, and on the order in as{fi), up to which the potential itself is expanded to, 
i.e. LO, NLO, NNLO and NNNLO. At NNNLO the couphng as{fi) in terms proportional 
to (//) is identified with ai^'°°''(/i), where I = 5 — n, and the Q!s(/i) in the logarithm 
with al~^'^°^{fi); at NNLO the identification works analogously, but Z = 4 — n, at NLO, 
I = 3 — n and at LO, / = 2 — n. The resulting expressions are functions of r and Aj^//Li. 
Consequently, Aj^ can be determined by a fit. 

Let us remark that both strategies (A) and (B) seem to be equally justified. There is no 
compelling reason favoring any of them. In particular both strategies ensure that at least all 
terms up to a given order in as(/u) are consistently taken into account. 

We perform these fits to various lattice results for the static potential, corresponding to the 
ensembles of gauge link configurations listed in Tabled! The fit is an uncorrelated minimizing 
two parameter fit of the perturbative expansion of the static potential V{r) to the corresponding 
lattice results in the range rmm 1^ r < rmax, where r^m, J'max and fi are input parameters, 
discussed below. The fit parameters are aVo (a constant shift of the static potential) and as(/i) 
(fitting procedure (A)) or oAj^ (fitting procedure (B)), respectively. 

4.2 Some general comments on our approach 

In passing, we shortly provide some general comments on our approach. This is intended to ease 
the reader in relating it to other works. 



12 



The first point to note, is that our approach is fully compatible with arguments emphasizing the 
necessity to account for renormalon contributions, in order to enhance the convergence behavior 
of the perturbative expansion. As outlined in detail in |38j . the expansion in terms of as(/i), 
evaluated at the fixed momentum scale /i, and the allowance for a constant shift of the static 
potential as independent fitting parameter at each order in the expansion of the QQ static 
potential, as employed by us, take the renormalon effects fully into account |38j . 

Second, a lot of work has been invested to determine the logarithmic contributions to the static 
potential. Those at ©(a^) have been determined explicitly in [33] • Moreover, resorting to 
the renormalization group, the leading ultrasoft logarithms [39] and next-to-leading ultrasoft 
logarithms [iQ] in the static potential could be resummed. This can be particularly important, 
when aiming at a study of the short-distance behavior of the static potential. In our work, we 
do not invoke the resummation of logarithms. Aiming at a direct comparison and fitting of the 
perturbative expressions with data points from lattice QCD simulations, we are naturally in the 
large-distance sector of the perturbative regime, where is not particularly small, but rather of 
the order of a few tenths. Consequently In Ug is of order 1 and, hence, not expected to dominate 
the other (unknown) contributions at a given order in the expansion in as- 

4.3 Input parameters to the fitting procedure and systematic errors of 
arising from the perturbative side 

4.3.1 Individual variation of input parameters 

In this section, we investigate the stability of A^ with respect to the input parameters rmin, 
''max and the scale /i, at which the coupling constant as(/u) is defined, using our finest lattice 
spacing (ensemble with (3 = 4.35). To this end we vary the input parameters in the following 
intervals: 

• ''min = 2a . . . 4a: 

lattice discretization errors in the static potential typically become small, when increasing 
the QQ separation beyond 2a . . . 3a. Due to the fact that our gauge action is tree-level 
Symanzik improved, and the lattice static potential is further improved, as explained in 
section 12. 2^ we expect only marginal discretization errors for r > 2a. 

• ''max = 4a. . . 6a (corresponding to rmax = 0.17fm . . . 0.25 fm): 

the NNNLO expression for the static potential is expected to be in good agreement with 
lattice data up to separations of around 0.25 fm [381 140j . 

• l/yLt = 3a . . . 5a: 

following the discussion in section 13.3.21 it is natural to choose l/^u in the vicinity of the 
fitting interval rmin ■ ■ -f max 

(cf. also [38]). 

Below, we demonstrate that the fit results for A^ are rather stable with respect to these param- 
eter variations, i.e. we confirm that a meaningful and a rather precise matching of perturbation 
theory and lattice results seems to be possible. 

Exemplary fits of the perturbative orders, i.e. LO, NLO, NNLO and NNNLO, to lattice results 
are shown in Figure [2j Fitting procedure (A) has been used. The above mentioned input 
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parameters have been chosen at the centers of then' above defined ranges of variation, i.e. r^i^ = 
3a, Tmax = 5a and = (rmin + ^max)/2 = 4a. The resulting values for A-^ values in MeV are 
shown in the figure captions. 
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P = 4.35 
1/fi = 4a, 



Exemplary fits of the perturbative expressions at LO, NLO, NNLO and NNNLO to 
lattice results. We employ fitting procedure (A), and set r^ain = 3a, rmax = 5a, and 



To understand, how the extracted depends on the input parameters rmin, rmax and 1//Li, 
we vary them individually in the following. 

In Figure [3] we vary the scale l//i between rmm = 3a and r^^x = 5a. Moreover, we compare 
fitting procedure (A) and (B): 

• In order to allow for a meaningful determination of it is important that the dependence 
of the fit results for Aj^ on the scale /i is rather weak; this amounts to a plateaux-like 
behavior in the plot of Aj^ as a function of l//i. Such behavior is clearly observable at 
NNNLO, and to some extent also at NNLO. 

• The expressions at NNNLO and NNLO yield similar results, indicating that higher orders 
do not contribute significantly. 
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The convergence of the four available orders, i.e. LO, NLO, NNLO and NNNLO, is better 
for fitting procedure (A) than for fitting procedure (B); this might be related to the fact 
that for (A) always the 4- loop expression for o^Cm) is used, i.e. that the fit formulae used 
for (B) contain less information. 

Fitting procedure (A): 

- Variation of A^jg (NNLO): 307 MeV . . . 366 MeV. 

- Variation of Aj^ (NNNLO): 312 MeV . . . 344 MeV. 

Fitting procedure (B): 

- Variation of Aj^ (NNLO): 257 MeV . . . 320 MeV. 

- Variation of Aj^ (NNNLO): 304 MeV . . . 313 MeV. 
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Figure 3: Plot of Aj^ as a function of l//i at LO, NLO, NNLO, and NNNLO; left: fitting 
procedure (A); right: fitting procedure (B). 



In Figure S] we investigate the stability of the fit results for Aj;^ with respect to a variation 
of the fitting interval rmin • • • J^max, applying fitting procedure (A); we vary both rmin and rmax 
separately as well as the center of the fitting interval, while keeping its extension r^ax — '"min 
fixed: 



• Fit results for depend only weakly on the fitting range in particular for NNLO and 
NNNLO. 

• Varying rmin = 2a ... 4a, r^ax = 5a, 1/fi = (rmin + ?'max)/2: 

- Variation of A,^ (NNLO): 327 MeV . . . 337 MeV. 

- Variation of Aj^ (NNNLO): 313 MeV . . . 332 MeV. 

• Varying rmax = 4a ... 6a, rmin = 3a, l/^a = (rmin + ?"max)/2: 

- Variation of Aj^ (NNLO): 325 MeV . . . 338 MeV. 
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- Variation of Aj^ (NNNLO): 316 MeV . . . 329 MeV. 
• Varying l//i = 3a . . . 5a, rmm = 1/^ - a, ^'max = 1/^ + a: 

- Variation of A,^ (NNLO): 318 MeV . . . 342 MeV. 

- Variation of Aj^ (NNNLO): 305 MeV . . . 336 MeV. 
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Figure 4: Aj^ at LO, NLO, NNLO, and NNNLO, fitting procedure (A); left: as a function of 
rmin; right: as a function of rmax! bottom: as a function of l/^u. 



The variation of Aj^ induced when varying the input parameters individually, as done in Figure [3] 
and Figure m is summarized graphically in Figure [5] both at NNLO (blue) and at (red). As 
systematic uncertainty one could quote, e.g. the whole range of values covered by the NNNLO 
variations, 

Aj^ = 304 MeV... 344 MeV. (46) 

The sources of the systematic error investigated above might, however, be correlated. To be 
on the safe side one could add the errors in quadrature. Another more sophisticated method is 
discussed in section 11.3. 21 
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Figure 5: Graphical summary of the variation of Aj^, when altering the input parameters 
individually (cf. Figure [3] and Figure [Hand the corresponding text). 

4.3.2 Inclusion of correlations between different sources of the systematic error 

To account for possible correlations, we perform a large number of fits, with the input parameters 
chosen randomly and uniformly in the intervals specified above. As systematic error we then 
take the variance of the fit results. 

We have performed 40,000 fits (sufficiently many that the statistical error of the variance is 
negligible) : 

• 10,000 NNLO fits, fitting procedure (A); 

• 10,000 NNLO fits, fitting procedure (B); 

• 10,000 NNNLO fits, fitting procedure (A); 

• 10,000 NNNLO fits, fitting procedure (B). 

For these fits we randomly choose 

• ^'min = 2a ... 4a and rmax = 4a . . . 6a, imposing the constraint rmax — ?^min > o,; 

• 1//^* ~ ^min • • • '"max' 

We obtain an average and a variance, i.e. a systematic error, of 

Ajjg = 315(26) MeV. (47) 

The error is slightly larger as compared to that stated in eq. (j46p . However, both results perfectly 
agree within errors. 
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The fitting procedure introduces an additional statistical uncertainty of ~ 2MeV, which is 
negligible, when added in quadrature. 

Further systematic uncertainties, exclusively originating from the lattice static potential, are 
discussed in the next subsection. 



4.4 Systematic errors of associated with the lattice computation 

In the following we address lattice discretization errors and finite volume effects. Moreover, we 
investigate the dependence of the static potential on the light quark mass. 

4.4.1 Lattice discretization errors 

To get an idea about the order of magnitude of lattice discretization errors, we turn to the lightest 
quark mass available at each of the four values of the lattice spacing. As can be seen from Tabledl 
the corresponding pions are roughly of the same mass ~ 284 MeV . . . 352 MeV. Keeping the 
parameters rmin and rmax in physical units approximately the same, we then extract values for 
Aj;^ in physical units, by performing fits (fitting procedure (A), NNNLO, l//i = (?'max+^min)/2), 

• /3 = 3.90: 

rmin = 2.00 X a, r^ax = 3.50 x a Aj^ = 317(12) MeV. 

• (3 = 4.05: 

ri^in = 2.53 X a, rmax = 4.42 x a A-^ = 312(10) MeV. 

• /3 = 4.20 {L^ xT = 48^ x 96): 

rmin = 3.14 X a, rmax = 5.50 X a A^j^g = 306(5) MeV. 

• /3 = 4.35: 

rmin = 3.77 X a, rmax = 6.60 x a Aj^ = 304(12) MeV. 

Note, that the errors associated with the extracted A^ values are rather large, because this time 
we have included the errors associated with the lattice spacings (cf. Table [T]). This is essential, 
when comparing results obtained at different values of the lattice spacing. 

These results for Aj^, as well as a continuum extrapolation assuming a dependence oc a^, 
are depicted in Figure EJ Within statistical errors, the result of the extrapolation, Aj^ = 
298(11) MeV, is in agreement with the result at /3 = 4.35. Nevertheless, Figure [6] indicates that 
there might be a slight downward tendency, when approaching the continuum limit. We are 
conservative and add an additional systematic error of ±6 MeV to Aj^, estimated by taking the 
difference between the central values of Aj^ at our smallest lattice spacing and the continuum 
extrapolation. 



4.4.2 Finite volume effects 

We investigate finite volume effects by considering the two spacetime volumes at /? = 4.20, which 
differ by a factor of 16. Again we extract A-^ by performing fits with identical parameters 
'"min = 3.14 X a and rmax = 5.50 x a: 
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Figure 6: Dependence of Aj-jg on the lattice spacing 



• (3 = 4.20 (L^ X T = 243 X 48): Ay^ = 305(1) MeV. 

• /3 = 4.20 {L^ xT = 48^ x 96): = 306(1) MeV. 

Within tiny statistical errors (compared to the systematic uncertainty of ±26 MeV determined in 
section 331) ^Ms does not change, when increasing the spacetime volume by the above mentioned 
factor 16. The smaller of the two lattice volumes investigated roughly corresponds to our 
ensemble at /? = 4.35, which we will use to determine the final result for A^. Therefore, 
we do not expect considerable finite volume corrections. 



4.4.3 Non-vanishing light quark mass 

Since the perturbative expressions in section [3] are derived under the assumption of massless 
dynamical quarks, while our lattice results employ light quarks of equal, but finite mass, with 
pion masses mps>284MeV, it is essential to investigate the quark mass-dependence. To this 
end, we consider the three different values of the light quark mass, available at /3 = 4.05 (cf. 
Tabled]). As before we extract Aj^ by performing fits with identical parameters rmin = 2.53 x a 
and Tmax = 4.42 x a: 

• mps = 325MeV: A^jg = 312(1) MeV. 

• mps = 449 MeV: A^ = 313(1) MeV. 

• mps = 517MeV: A^jg = 312(1) MeV. 

Within tiny statistical errors Aj^ is constant in the quark mass region investigated. Therefore, 
we do not expect a dramatic change, when approaching the massless limit. In other words we 
consider the systematic error introduced by comparing a massive lattice computation with a 
massless perturbative calculation negligible compared to the uncertainty of ±26 MeV already 
determined in section 14.31 
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4.5 Final results for 

In the following, we quote our final results for in QCD with nj = 2 dynamical quark flavors. 
We list these results both in units of rc0 and in MeV. The determination is based on the lattice 
results at our finest lattice spacing, i.e. at /3 = 4.35, as explained in section 14.3.21 The errors, 
which we combine by adding them in quadrature, are 

(1) the systematic errors associated with perturbation theory and the input parameters of the 
fitting procedure (cf. section [4.3.2p . 

(2) the estimated lattice discretization errors (cf. section [4. 4. Ih . 

(3) the errors associated with r^/a = 9.81(13) and the lattice spacing a = 0.0420(17) fm, 
respectively. 

Therewith, we finally obtain 

roA^jg = 0.658(55) , = 315(30) MeV. (48) 



''ro is defined via roF{ro) — 1.65, where F{r) — dV(r)/dr gT]. 
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5 Conclusions 



We have determined Aj^jg for QCD with nf = 2 dynamical quark flavors, by fitting various orders 
of the perturbative expansion of the QQ static potential, up to terms of ©(af) and ~ Ina^, 
to corresponding lattice data points. 

In units of the hadronic scale rg our result reads 

roAj^ = 0.658(55), (49) 
while in physical units it is given by 

Aj^ = 315(30) MeV. (50) 

To obtain the latter result, the physical scale has been set by fixing the lattice scale via the 
experimental results for the pion mass and the pion decay constant 



f 



All sources of systematic error have been investigated: neglect of higher orders in the perturba- 
tive expansion, the choice of the scale fi, at which Os is defined, the dependence of the fit results 
on the fitting range, lattice discretization errors and finite volume effects, as well as finite quark 
masses. The errors quoted above are dominated by the variations of Aj^, when fitting different 
perturbative expressions for the static potential to the lattice results, while varying the fitting 
window and the scale fi. In principle this error could be reduced. This, however, would require 
either higher orders in the perturbative expansion, or lattice results at even finer lattice spacing; 
neither is available at the moment. 

Our results compare well with other determinations of Ay^ for QCD with nj = 2 dynamical 
quark flavors from the literature, e.g. 

• a lattice computation using the Schrodinger functional [131111]; result: ?^o^ms ~ 0.73(3) (5); 

• a lattice computation using tq and a boosted coupling [l5]; result: Aj^ = 261 (17) (26) MeV 

(the scale has been set via = 0.467 fm, while our ETMC rp is significantly smaller, 
j^ETMC ^ 0.42 fm; converting the result to our scale yields A^jg = 290(19)(29) MeV); 

• a lattice computation using Landau-gauge gluon and ghost correlations [l6]; result: tqA^ = 
0.60(3)(2); 

• an ETMC lattice computation using the ghost-gluon running QCD coupling [U]; result: 
Aj^ = 330(23) (22)_33 MeV; 

• variationally optimized perturbation, combined with renormalization group properties |48j : 
result: Aj^ = 255+^^ MeV. 

Our result for A^ and the results from literature listed above are summarized graphically in 
Figure El 



^Since the lattice setup is n/ = 2 QCD, while experimental results may also be affected, e.g. by dynamical 
strange and charm quarks, isospin breaking or electromagnetic effects, a different choice of observables, for example 
the pion mass and the nucleon mass, might yield a slightly different lattice spacing and, therefore, a corresponding 
change of the result for in units of MeV. The effect might be as large as 10%, i.e. of the same order of 

magnitude as the error ±30 MeV already quoted. Therefore, specifying in units of ro is preferable. For a 
more detailed discussion about problems and ambiguities associated with scale setting cf., e.g. [42| . 
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comparison of various Attc results 
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Figure 7: Graphical summary of various results for in QCD with nj = 2 dynamical quark 
flavors (to eliminate discrepancies due to different methods of scale setting, the results from 
|43 1 144 ^ HSj I46j have been converted using our value for rg, Tq^^*-' = 0.42 fm). 



Let us finally emphasize again that our value for for QCD with n/ = 2 is in good agreement 
with the results of other approaches aiming at its determination (cf. Figure [7]). The systematic 
error in our final result, eq. (j50p . is of the order of ten percent, which is also compatible with 
the systematic errors of the works listed in Figure [71 Our investigation here opens the path 
towards the determination of Aj^ for a larger number of dynamical quark flavors also, such 
as the ongoing efforts of the ETMC collaboration using Nj = 2 + 1 + 1 (dynamical up and 
down, strange and charm quarks) [491 150|. As we were able to show in this work, that Aj^ can 
be determined from the static potential for a lattice spacing of a ~ 0.042 fm, a corresponding 
calculation for Nf = 2 + 1 + 1 will hence be possible as soon as respective lattice data at such 
small values of the lattice spacing becomes available. 
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A Monte Carlo evaluation of the integral ([9]) 



In order to achieve a tree-level improvement of the static potential, computed with the tree-level 
improved Symanzik action, one has to solve the integral 

air) = -I-TA— , nl^icosfe^,) — — 

(cf. section [2] and eq. Q). 

In principle one could evaluate this integral directly by means of standard Monte Carlo methods. 
However, due to the fact that the integral is dominated by the singularity at k = 0, this would 
result in a rather large statistical error. 

To circumvent this problem, we split the integral in two parts, 

G(r) = Gi(r) + G2(r) (52) 
1 r%3. n?=iCos(A;,r,) 



""'^'^ - i^L ^^iC^ '''' 

1 r\3r f n?=iCOs(%rj) lfj=iCOsikjr. 



(27r)3y_^ " ■A4E|=iSin2(A;,/2) + (4/3)E|=iSm^(A:,/2) 4 E?=i sin2(A;,/2) 

(54) 



The first part, Gi, can be solved exactly by applying a clever recursion technique [51j . Since the 
integrand of the second part is rather smooth, a Monte Carlo evaluation yields a statistical error, 
which is around two orders of magnitude smaller than that obtained by a direct computation 
of ()5ip . when using a comparable amount of computer time. In other words the necessary 
computational resources needed to reach a desired precision are reduced by a factor of « lO'^. 
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